Phase separation of crystal surfaces: a lattice gas approach 
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We consider both equilibrium and kinetic aspects of the phase separation ("thermal faceting") 
of thermodynamically unstable crystal surfaces into a hill-valley structure. The model we study 
is an Ising lattice gas for a simple cubic crystal with nearest-neighbor attractive interactions and 
weak next -nearest-neighbor repulsive interactions. It is likely applicable to alkali halides with the 
sodium chloride structure. Emphasis is placed on the fact that the equilibrium crystal shape can 
be interpreted as a phase diagram and that the details of its structure tell us into which surface 
orientations an unstable surface will decompose. We find that, depending on the temperature and 
growth conditions, a number of interesting behaviors are expected. For a crystal in equilibrium with 
its vapor, these include a low temperature regime with logarithmically-slow separation into three 
symmetrically-equivalent facets, and a higher temperature regime where separation proceeds as a 
power law in time into an entire one-parameter family of surface orientations. For a crystal slightly 
out of equilibrium with its vapor (slow crystal growth or etching), power-law growth should be the 
rule at late enough times. However, in the low temperature regime, the rate of separation rapidly 
decreases as the chemical potential difference between crystal and vapor phases goes to zero. 

PACS Numbers: 68.35.Rh, 68.35.Bs, 64.60.Cn, 64.60.Ht 
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I. INTRODUCTION 

The decomposition (or "faceting") of a surface into 
pieces of surfaces of other orientations is a problem 
of long-standing interest in the materials science and 
physics communities. It can occur in a variety of cir- 
cumstances, including when a crystal is being grown or 
etched [Q, when deposition is occurring through such 
processes as molecular beam epitaxy [^|,||, when other 
materials are adsorbed onto the surface H, or when an 
electric current is applied across the sample How- 
ever, such external perturbations are not always neces- 
sary to cause this decomposition: If the surface tension is 
sufficiently anisotropic, some surface orientations will be 
thermodynamically unstable |Q. Then, a surface which 
is initially prepared in such an orientation will spon- 
taneously decompose into a faceted structure, evolving 
towards a final equilibrium configuration consisting of 
large facets of stable orientations. This decomposition 
has been variously referred to as "equilibrium faceting," 
"thermal faceting," "thermal etching," "Herring recon- 
struction," or "hill- valley reconstruction" 

Recently, there has been a resurgence of interest in this 
process in both the experimental ||l(]|-jl2fl and theoretical 
physics communities [|l3|-[T^]. In the theoretical commu- 
nity, the close analogy to phase separation of binary liq- 
uids or alloys has been pursued. Hence, the process has 
been dubbed "phase separation" or "spinodal decompo- 
sition" of crystal surfaces 14 1^]. (The terms "faceting," 
"phase separation," "decomposition," and "coarsening" 
each have their drawbacks in describing this process in 
all its various manifestations but, lacking better termi- 
nology, we will use all these terms and use them more or 
less interchangeably.) 



The process, which is illustrated in Fig. [|, proceeds 
as follows: The surface, initially prepared in an unstable 
orientation, first decomposes into small pieces ( "facets" ) 
of stable surface orientations, in spite of the fact that 
this increases the total surface area. The faceted surface 
then coarsens over time in order to minimize the energy 
associated with the edges between the different surface 
orientations. The dynamical process by which the size of 
the facets grows over time is indeed closely analogous to 
the phase separation of a binary liquid or alloy, and recent 
theoretical work has attempted to elucidate the analogy 
and to investigate the dynamics associated with such sep- 
aration p3|-p6[| . Much earlier work along these lines was 
performed by Mullins ||, and he predicted power law 
growth of the facet size L with time t: L(t) ~ t™ with 
n = 1/4, 1/3, or 1/2 for the mass-transfer mechanisms 
of surface-diffusion, volume-diffusion, or evaporation- 
condensation, respectively. However, the experimental 
situation has remained murky. Generally, it has been 
found that such faceting occurs only very slowly {e.g., 
n « 0.1) if at all under near-equilibrium conditions (see 
discussion in |ll|), and indeed, some have argued that 
such faceting is not even a thermodynamic phenomenon 
but occurs only under driven conditions of crystal growth 
or etching [[[J. Nonetheless, it is known theoretically that 
such faceting should in principle occur in equilibrium if 
the surface tension is sufficiently anisotropic and some 
of the experimental observations do support this 
Thus we are led to conclude that one should seek a dy- 
namical explanation to understand why so little faceting 
is seen 
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Two of the recent studies of the dynamics have used 
solid-on-solid approximations of an Ising lattice gas on 
a cubic lattice with interactions appropriate for alkali 
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halide materials such as sodium chloride (NaCl) Jl9| . [20j . 
Shore, Holzer, and Sethna (b|] considered the special 
case of a [111] surface and argued that the coarsening 
of the surface should be only logarithmic with time (at 
long times) for quenches to all temperatures T at which 
the [111] surface is unstable. This claim was supported 
by their Monte Carlo simulations. Vlachos, Schmidt, 
and Aris |l6| performed Monte Carlo simulations for the 
coarsening of [hkO] surfaces, concentrating in particular 
on a [210] surface. Although they draw no firm conclu- 
sions about the asymptotic growth laws, we believe that 
the results presented there — of anomalously low expo- 
nents which decrease with decreasing temperature — are 
compatible with (although certainly not proof of) loga- 
rithmic growth at long times [ |13|]2l|l . 

The approach of Stewart and Goldenfeld jl4| and Liu 
and Metiu [jl5| has been along somewhat different lines, 
closer in spirit to the original work of Mullins. They have 
looked at the problem in the continuum limit, which is 
generally claimed to be valid when one is interested in 
behavior at long times and on long length scales. In this 
limit, the driving force for the initial breakup of the sur- 
face into facets is provided by a sufficiently anisotropic 
surface tension. Also included in the free energy func- 
tional is a term which suppresses rapid changes in the 
surface normal |Ts|| . Such a term has two effects: First 
it rounds off the edges and corners on the surface, thus 
preventing the occurrence of singularities, which would 
lead to divergences in the Langevin equations. Second, 
it reproduces, in a macroscopic sense, the energy costs 
associated with edges and corners, thus providing the 
driving force for the coarsening of the surface structure. 
For this model, Liu and Metiu concluded that coars- 
ening proceeds logarithmically for quasi-one-dimensional 
ordering. Simulations of the resulting Langevin equation 
for a two-dimensional surface gave L{t) ~ t n with n w 
0.13 and 0.23 for surface diffusion and for evaporation- 
condensation mechanisms, respectively. These results 
were noted to be in fairly good agreement with sim- 
ple power-counting arguments on the Langevin equation, 
which suggest n = 1/6 and 1 /4 for the two different mech- 
anisms. 

Although continuum approaches are a powerful tool, 
and currently enjoy great favor in the study of kinetics of 
growth, they suffer from at least two major limitations 
that make it useful to also consider microscopic {e.g., 
lattice gas) models. The first is that, since they deal with 
the system on a coarse-grained level, one generally needs 
to assume some form for the input parameters, in this 
case the orientation-dependence of the surface tension. 
Thus, microscopic models are useful to give us guidance 
on what form to choose for the surface tension. Indeed, in 



Sec. Ill B we will find that the form of the surface tension 



that arises in one such microscopic model, and may occur 
quite generally is considerably different from that which 
one might have expected from naive considerations. 

The second limitation is that such models, at least as 
currently formulated, appear to be unable to properly 



model aspects of the problem where the atomic discrete- 
ness is fundamental and cannot be ignored. This leads 
to the conflicting predictions discussed above (logarith- 
mic vs. power-law growth) between the microscopic and 
continuum approaches for the coarsening of an unstable 
[111] surface of NaCl. 

The lattice gas models, of course, have their own lim- 
itations {e.g., the neglect of elastic effects Jl^]), but at 
this point a fairly complete picture of the equilibrium 
and kinetic aspects of the faceting process within a mi- 
croscopic model would be useful in advancing the under- 
standing of the faceting process in the real world. In this 
paper, we will try to present such a unified picture of 
faceting within the context of the Ising lattice gas model 
for NaCl. The picture combines and generalizes some of 
our earlier work on the faceting of a [111] surface [ p"3| , 
and some more recent results concerning the phase di- 
agram of the six-vertex model p3|. We also include a 
discussion of how the kinetics of the faceting process will 
be altered for the case of a driven interface, that is, one 
in which the crystal is being slowly grown or etched. 

The outline of the paper is as follows: In Section II, 
we introduce the Ising lattice gas model which we will 
study, discuss the general features of its equilibrium crys- 
tal shape (ECS), and explain how the ECS tells us which 
surface orientations are unstable and into what surfaces 
such an unstable orientation will decompose. In Section 
III, we consider in detail the equilibrium and kinetic as- 
pects of the phase separation in two different tempera- 
ture regimes. In Section IV, we discuss how the kinetics 
is altered in the case when the crystal is no longer in 
equilibrium with its vapor (growth or etching). Finally, 
in Section V, we present a summary of our conclusions. 



II. EQUILIBRIUM CRYSTAL SHAPE 

The problem of determining the equilibrium shape of 
a crystal by minimizing its total surface free energy sub- 
ject to a constraint of fixed volume has been considered 
for over a century p3| . The major result is the cele- 
brated Wulff construction by which one determines the 
ECS given the surface free energy (surface tension) as a 
function of the surface orientation. 

Only much more recently, however, has the ECS of 
some simple microscopic models for crystals been consid- 
ered p9lp4-26|. Of particular interest to us is the Ising 



lattice gas on a simple cubic lattice with nearest-neighbor 
attractive interactions and next-nearest neighbor repul- 
sive interactions, first studied by Rottman and Wortis 
p9|| . The Hamiltonian is 
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where the Si take on the values +1 and —1. In the lat- 
tice gas context, +1 is taken to represent an atom oc- 
cupying the site and —1 is taken to represent an empty 
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site. The first sum is over all pairs of nearest-neighbors 
(NN) while the second is over all pairs of next-nearest- 
neighbors (NNN). We have chosen our sign convention 
so that both J\ and J2 are positive when the NN bonds 
are attractive and the NNN bonds are repulsive. We 
are interested in the case where J1/J2 > 4, in which case 
the ground state is (using spin language) "ferromagnetic" 
& 

In many instances, it will be useful to consider an 
interface between occupied and unoccupied sites (i.e., 
between solid and gas) in the_[lll]-restricted-solid-on- 
solid (RSOS) approximation p0||l3j ]. This approximation 
can be obtained by viewing the interface from the [111] 
direction and requiring that none of the interface be hid- 
den from view by other parts of the interface (see Fig. [l]). 
It is equivalent to taking the limit J\/ J2 — ► 00 with T / J 2 
fixed. Configurations in this model can also be viewed 
as tilings of the plane by 60° rhombi of three different 
orientations. In this representation, the energetics are 
reproduced by assigning an energy of 2 J2 to each border 
between unlike rhombi, and the coarsening process in- 
volves the phase separation of the three types of rhombi. 

Finally, A/z = [x s — fj, v represents a difference in chemi- 
cal potential between solid and vapor. (In spin language 
it would be a magnetic field.) For the remainder of this 



paper, with the exception of Section IV, we will be con- 
sidering the case in which the solid and vapor are in equi- 
librium (A/i = 0). 

The model of Eq. (|l|) and its RSOS approximation were 
proposed to represent materials such as sodium 

chloride (NaCl), where the ions of different species (and 
opposite charge) would have an attractive interaction 
while those of the same species would repel one another 
|3. The ECS for this model is show in Fig. § It 
was found Q that the ECS remains strictly cubical, 
with (macroscopically) sharp edges and corners up to the 
corner-rounding transition at a temperature T = Tcr, at 
which point the crystal first rounds at the corners. For 
the RSOS model, Tqr can be calculated exactly [ ^0|Jl3| 
and is given by pq] 



T, 



-4J 2 



CR 



ln[l/3-5/(9a 1 /3) +Q! i/3] 



7.1124... J 2 



(2) 



where a = + \/23/3). As the temperature is further 
increased, the rounded region spreads out along the edges 
until the edge-rounding temperature Ter, at which point 
the entire edge becomes rounded. The smooth {100} 
facets still remain, however, up to the roughening tem- 
perature Tr p9[] . Finally, above Tr, the entire crystal 
shape is rounded. 

We make two further observations about these inter- 
facial phase transitions. The first is that when J2 = 
or is attractive, then Tcr = Ter = 0. Therefore the 
corner and edge rounding transitions are the result of 
the competing attractive and repulsive interactions. The 
second is that within the [111J-RSOS model, Ter and 
Tr are infinite. This is a result of the fact that setting 



Ji — > 00 suppresses the fluctuations which are responsi- 
ble for edge rounding and the roughening of the facets. 
Thus, the usefulness of the RSOS approximation is re- 
stricted to temperatures T < Ter. It is most reliable for 
surface orientations close to [111]. 

Experimentally, in NaCl, the corner-rounding transi- 
tion is found to occur at a temperature Tcr ~ 920 K 
|]l7| , p0 24, 3(]]. The sharp cubical shape and the shape 
with facets but no sharp corners or edges have both been 
observed experimentally [T^] . The intermediate shape 
with sharp edges but rounding at the corners was also 
seen on some crystals. However, it is likely that these 
crystals were not sufficiently equilibrated, so one would 
have to say that there is, as of yet, no experimental con- 
firmation of this intermediate shape, and thus, of the 
claim that the corner- and edge-rounding temperatures 
are distinct (See Jl7| and the discussion in [[24f). 

It is a theorem due to Herring that those surface 
orientations which do not appear on the ECS are thcrmo- 
dynamically unstable (and will thus break up into pieces 
of surface of stable orientations) . If we represent each sur- 
face orientation by a unit vector normal to that surface, 
then we can use a unit sphere to represent the various 
surface orientations. For the model specified by Eq. (ffl), 
Fig. |H shows an eighth of such a unit sphere, with shading 
on the sphere used to indicate those surface orientations 
which are stable in various temperature regimes. Since 
only the {100} facet orientations appear on the crystal 
shape (Fig. ^) for T < Tcr, they are the only stable 
surface orientations in this temperature regime. At Tcr, 
the [111] surface orientation (which is at the center of 
our diagrams in Fig. |§|) and all those orientations on the 
lines connecting it to the three nearest facet orientations 
become (marginally) stable. As the temperature is raised 
further, more and more surface orientations become sta- 
ble. However, it is not until the temperature reaches Ter 
that all surface orientations appear on the ECS and are 
thus stable. 

Herring's theorem is actually a specific consequence of 
a more general fact: The ECS is essentially a phase di- 
agram for the surface orientations 31 24|^|. More pre- 
cisely, the height of the ECS above a reference plane, 
F(h,v), as a function of the two in-plane coordinates 
(h, v), can also be considered as a free energy as a func- 
tion of "fields" (h,v) p2|. These fields are conjugate to 
the two-component order parameter (a;, y) specifying the 
surface orientation. The projection of the surface tension 
onto the reference plane is then given by the free energy 
as a function of surface orientation F(x, y). As is usually 
the case, the conjugacy of the fields and order parame- 
ters means that these two free energies are related by a 
(two-dimensional) Legendre transformation [ pl| : 



F(h,v) — min{T(x,y) — hx 



(3) 



which implies F(h, v) = F(x, y) — hx — vy with x and 
y chosen such that dF(x,y)/dx = h and dF(x,y)/dy = 
v. Equivalently, we can write x = —dF(h,v)/dh and 
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y = —dF(h, v)/dv, which simply states the fact that the 
surface orientation is specified by the slope on the crystal 
shape F(h, v). 

Because of the mapping between the ECS and a free 
energy surface, we can transfer all our terminology and 
knowledge about phase diagrams over to the realm of 
equilibrium crystal shapes! For example, sharp edges 
in the ECS are places where F(h, v) has a discontinu- 
ity in the first derivatives, i.e., first-order phase transi- 
tions; whereas boundaries where a facet smoothly joins 
the rough regions of the crystal are places where the 
first derivatives of F(h,v) are continuous but higher- 
order derivatives are not, and they are thus considered 
second-order transitions. Even the scaling behavior at 
these boundaries can be analyzed using the tools of crit- 
ical phenomena p6p^ ]. 

At first-order lines ("edges"), one has coexistence be- 
tween the two "phases" (surface orientations). The cor- 
ners of the ECS for T < Tcr are points of coexistence 
between three orientations. Such an interpretation of 
edges and corners as coexistence lines and points is im- 
portant because it means that once we know the details 
of the ECS, we can immediately determine into which 
surface orientations an unstable surface will decompose. 
The consequences of this will be further spelled out in 
the next section where we consider the phase separation 
process in detail in the two different temperature regimes 
T < Tcr and Tor <T < T E r. 

III. PHASE SEPARATION FOR THE CASE OF 
SOLID-VAPOR EQUILIBRIUM 

We first look at a few snapshots of equilibrium config- 
urations of a surface of orientation [15, 5, 1] in the [111]— 
RSOS model. In Fig. |(a), we see that for T < T CR the 
configuration consists of a surface which (apart from a 
few thermal fluctuations) is phase-separated into [100], 
[010], and [001] facets. In Fig. ^(b), at a temperature 
T CR < T < T ps ([15,5,l]) where T ps ([15, 5, 1]) is the 
phase-separation temperature for this surface, it appears 
to have separated into one phase vicinal to the [100] facet 
and another vicinal to the [010] facet. Finally, at a tem- 
perature T > T ps ([15, 5, 1]), the surface is stable and does 
not phase separate. These qualitative observations set 
the stage for our detailed study of the equilibrium and 
dynamic aspects of the phase separation problem in the 
two regimes, T < T CR and T CR < T < T ps ([hkl}) < T ER , 
in which phase separation occurs for an arbitrary surface 
[hkl]. 

A. Below Tea 

In the regime T < Tcr, all surfaces except the {100} 
facets are unstable. Any surface prepared in another 
orientation will decompose into a combination of these 



surfaces. We can determine which surfaces by utilizing 
the mapping of the ECS onto a phase diagram discussed 
above. For example, an arbitrary surface [hkl] with h, k 
and / all positive is a surface which "lives" in the three- 
phase coexistence region represented by the corner of 
the ECS where the [100], [010], and [001] facets meet. 
Therefore, such a surface will decompose into a linear 
combination of these three surface orientations, with the 
amount of each orientation determined by the require- 
ment that the resulting surface still has the average ori- 
entation [hkl]. (Electron micrographs of the faceting of 
the [111] surface of NaCl can be found in Ref. ||, but see 
also the comments in Ref. JT7|].) For the special case of a 
[hkO] surface, the surface "lives" in the two-phase coex- 
istence region represented by the entire edge between the 
[100] and [010] facets and therefore decomposition occurs 
into just these two surface orientations. 

In Ref. |l3| , we gave arguments for logarithmically slow 
growth of the facet sizes below Tcr for the case of the de- 
composition of a [111] surface. There, we also presented 
strong numerical evidence from Monte Carlo simulations 
supporting this claim. Although the discussion there was 
for evaporation-condensation dynamics and within the 
[lll]-RSOS model, the same arguments should hold out- 
side the RSOS approximation and should also apply for 
the case when the dominant mechanism is surface diffu- 
sion (in which case the RSOS model is too restrictive to 
be used). The arguments also generalize to the case of 
decomposition of any arbitrary surface [hkl] . 

While the reader is referred to Ref. |13| for the numer- 
ical evidence supporting our claim, the basic argument 
itself is repeated here for completeness: Consider a coars- 
ening surface such as that shown in Fig. [I]. At a time t, 
the characteristic length scale (i.e., average facet size) 
is L(t). In order for the structure to coarsen further, a 
step across must propagate across one of the facets; for 
definiteness, say a [100] facet. Since the step consists 
of, microscopically, pieces of [010] and [001] surface ori- 
entations embedded in the [100] surface orientation, the 
step consists of (again, on the microscopic level) edges 
between facets of different orientation, and thus it costs 
an energy per unit length. Once the effects of entropy are 
considered, we find that there is still a nonzero step free 
energy (i.e., a free energy per unit length) up to T = Tcr. 
In fact, within the RSOS approximation, this step free 
energy as a function of angle can be calculated exactly 
(See appendix). It is only at Tcr that this step free en- 
ergy goes to zero for a step of a certain angle (specifically, 
a 45° step). 

The coarse-grained differential equation for the growth 
of L(t) during the coarsening process can be written as 

H 

dL = a(L, T) 

dt L m 1 ' 

Here, m depends on the specific type of dynamics to be 
considered. For example, for evaporation-condensation 
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dynamics, we likely have m = 2 or 3 (See Sec. B). How- 
ever, the result below for the asymptotic growth law is 
independent of m. 

In the standard case where the kinetic coefficient 
a(L,T) has no (or only very weak) dependence on L, 
Eq. d) gives L(t) ~ tV^+i), However, here, the fact 
that the coarsening involves activation over energy bar- 
riers of height fs(T)L implies that the kinetic coefficient 
a(I/,T) decreases exponentially with L: 



a(L,T) = a e-f^ L / T 



(5) 



/b(T) is a free energy barrier per unit length, and is 
related to the step free energy as discussed at the end 
of the appendix. At asymptotically long times, Eqs. (ji|) 
and (H) imply 



L(t) 



T 



fs(T) 



ln(t) 



(0) 



It is worthwhile to consider this result in the context of 
the recent work by other groups on this problem. Monte 
Carlo simulations of the coarsening process for [hkO] sur- 
faces, using essentially the same model as we have dis- 
cussed above, have recently been presented by Vlachos, 
Schmidt, and Aris p6| . They consider both the case of 
the crystal in equilibrium with the vapor and the case 
where there is a chemical potential difference between 
them. For the equilibrium case, they found generally slow 
coarsening with an effective exponent n e s in L(t) ~ t n ° a 
which is small and decreases with decreasing tempera- 
ture. This scenario is consistent with the behavior found 
in previous studies of models believed to coarsen loga- 
rithmically |L3|,^l|. Since small exponents are hard to 
distinguish from a logarithm (and the logarithmic form 
is expected only asymptotically), the ability of Vlachos 
et al. to fit their numerical data to power law growth 
with a small exponent should not surprise us. 

We should note that the claim by Vlachos et al. that 
a logarithmic form is inconsistent with their data fllf l 
is intended to mean only that a logarithm cannot be fit 
over the entire time regime f34| . In fact, from the inset 
of their Fig. 2(b), we see that a logarithm fits quite well 
(at least as well as a small power law) over the later 
time regime (t > 10 2 ). The data for shorter times are 
in the so-called "facet nucleation regime" |l6| and can 
be fitted by neither a logarithm nor a power law. This 
regime can be understood as the characteristic activation 
time for depositing or desorbing particles on the initial 
surface, since the smallest energy barrier of 4Ji implies 
a characteristic time on the order of t ~ exp(4Ji/T) ~ 
55 to successfully deposit the first layer on the surface. 
(Here we use the fact that the energy scale w\ in Ref. 
[ PI is related to ours by w± = 4Ji.) 

Other results presented by Vlachos et al., such as the 
speeding up of the growth rate when a chemical poten- 
tial difference exists between crystal and vapor, are at 
least in qualitative agreement with our picture (See Sec. 
IV). Thus, we can say that the simulations performed 



in Ref. |iq ] , while certainly not providing very conclusive 
support for our claims, are at least consistent with them. 

The continuum approach of Liu and Metiu p5| , by 
contrast, clearly predicts power-law growth of the facet 
sizes. How are we to reconcile this difference? Our as- 
sertion is that this continuum approach misses a funda- 
mental part of the problem which proves vital to the dy- 
namics in the discrete models (and, we believe, although 
we are somewhat less certain, in the experimental sys- 
tems). In particular, in those cases where an unstable 
surface decomposes into surfaces which are smooth (i.e.. 
below their roughening transition), the dynamics of the 
faceting of surfaces generically involves the propagation 
of a step across the facet. Since such a step has a nonzero 
free energy per unit length for T < Tcr, this leads 
to length-scale-dependent barriers and logarithmically- 
slow growth. 

Such a dependence is not captured by a continuum 
model which does not recognize that there is a smallest 
size, i.e., the width of a step, determined by the dis- 
creteness of the system. Rather, in Ref. |D| the authors 
introduce a term into their free energy which imposes an 
energy cost for rapid changes in surface orientation. This 
captures some of the physics that arises in the micro- 
scopic models, but clearly not all of it. Furthermore, be- 
cause the Langevin equations cannot cope with singular- 
ities (and in keeping with the spirit of coarse-graining), 
they round the cusps in the surface tension associated 
with the crystal facets. This means, in effect, that there 
is no longer any roughening transition in the model: all 
surface orientations are rough! 

The idea that the discreteness of a crystal is fundamen- 
tal in determining dynamical, and even equilibrium, be- 
havior is certainly not without precedent. The very exis- 
tence of a roughening transition and of nucleated dynam- 
ics for crystal growth below this transition is dependent 
upon including terms in the Hamiltonian (and the re- 
sulting dynamical equations) which explicitly model the 
discreteness of the system |35|-|3q| . As far as we know, 
there is no known prescription for determining a priori 
whether such discreteness will be relevant or irrelevant in 
a renormalization group sense. In the absence of such a 
prescription, purely continuum approaches to the prob- 
lem (even if they mimic some effects produced by the 
discreteness, e.g., by using anisotropic surface tensions 
and energy penalties for changes in surface orientation) 
should, we believe, be used with a certain degree of cau- 
tion. 

In Sec. IV, we will return to the case of coarsening be- 
low Tcr, but in the case where the crystal is slightly out 
of equilibrium with the vapor. There, we will find that 
the difference in chemical potential between vapor and 
crystal sets a maximum size to the free energy barrier, 
so the problem becomes one of nucleation and there is a 
return to power-law behavior at late enough times. How- 
ever, first we will consider the coarsening process for the 
case of crystal-vapor equilibrium but in the temperature 
regime above Tcr- 
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B. Above Tor 



alternately, along the sharp edge between the two facets) 



When we consider the diagram in Fig. || showing the 
stable orientations at a temperature Tor < T < Ter, 
three questions immediately present themselves: (1) How 
does one determine the location of the boundaries be- 
tween the stable and unstable surface orientations? (2) 
For a surface in the unstable region, what are the orien- 
tations of the stable surfaces it breaks up into? (3) What 
are the kinetics of this process? 



1. Equilibrium aspects 

We find it most natural to begin with a discussion of 
the second question raised above. The behavior of an 
unstable surface depends on the details of the ECS. In 
Figs. H and ||, we consider two possible scenarios 

In Fig. ||(a), for what we dub the "ridge scenario", the 
sharp edges ("ridges") between the facets continue into 
the curved region of the ECS. As one moves along such a 
ridge into that region, the angle between the unit vectors 
normal to the two surfaces that meet at that point on the 
ridge decreases from 90° at the facets to 0° at a second- 
order point at which the ridge terminates. Since at any 
point along this ridge two surface orientations meet, the 
mapping of the ECS onto a phase diagram tells us that 
this means there is coexistence between two surface ori- 
entations at each point. Without loss of generality, we 
consider the case of a surface [hkl] with h > k > I > 0. 
Fig. ||(b) , showing part of the diagram of stable surface 
orientations, demonstrates how such a surface prepared 
at an unstable orientation breaks up into two surfaces. 
Note that, by symmetry, the two coexisting surfaces are 
symmetric about the line h = k. Also, a linear com- 
bination of these two surface orientations must add up 
to give an average surface orientation of [hkl]. These 
two conditions, along with the third condition specifying 
those surfaces which are on the curve of marginally sta- 
ble surfaces, completely specify the two surfaces which 
the surface [hkl] decomposes into. 

In Fig. ^|(a), which we dub the "conical point scenario" , 
the sharp edge terminates at the point where it meets the 
facet boundaries. (This can be considered as a limiting 
case of taking the length of the ridges into the curved 
regions to zero.) In this case, all the unstable surfaces 
must "live" in the point where the edge meets the facet 
boundaries. Furthermore, all the marginally-stable ori- 
entations along the boundary in Fig. ^(b) must meet at 
this one point on the ECS. The slope of the ECS at 
this point depends on the angle from which it is entered. 
Therefore, in the rounded region of the ECS, the point 
has the symmetry of the tip of a cone, hence the name 
"conical point." Also note that the ECS has a jump in 
slope as one goes from the rounded region of the ECS, 
through the conical point, and onto one of the facets (or, 



Naively, one might expect that the ridge scenario 
would occur in general, as the conical point scenario 
seems to require that there be special symmetry about 
the conical point, which results in the meeting of all the 
different orientations at this sing le point on the ECS 
Such a point is very much analogous to the zero field 
(H = 0) point in the low temperature phase of a three- 
dimensional XY model, at which an entire one-parameter 
family (circle) of magnetizations coexist. However, in 
that model such degeneracy is clearly the result of the 
symmetry with respect to spin orientation 9 which oc- 
curs for H = in the original Hamiltonian. In the 
present problem, there is no such obvious symmetry and 
thus such a point would not be expected unless the terms 
breaking this symmetry become irrelevant under renor- 
malization of the original Hamiltonian. 

To rigorously determine which scenario occurs in one 
specific model, we have recently studied plj the [110]- 
RSOS model for the surface of an FCC crystal with 
nearest-neighbor attractive and next-nearest neighbor 
repulsive interactions p6|| . In this model, Tqr = 0, i.e., 
the corners round for any T > 0. However, a sharp 
edge separating the [111] and [llT] facets persists up 
to a nonzero temperature Ter. The advantage of this 
model is that it can be mapped onto the six-vertex model 
H,|H which is exactly solved |p3p^|. The solution of 
the six-vertex model, which is given in terms of integral 
equations, can be investigated numerically and, in cer- 
tain limits, analytically p^ ]. We find that in this model, 
the conical point scenario occurs, as was conjectured by 
Jayaprakash and Saam [ p6| at the time when they first 
considered this RSOS model f||. 

We also find an exact expression for the curves sepa- 
rating the stable and unstable surface orientations, as a 
function of temperature. (These curves are analogous to 
the boundaries between the shaded and unshaded regions 
for Tcr < T < Ter shown in Fig. || for the simple cubic 
case). As one moves along the ECS away from the con- 
ical point into the curved region, the crystal shape has 
the following nonanalytic behavior: 



F{h) = F + F x 5h + F 2 Sh 5/3 + 0{6h 7/3 ) 



(7) 



where Sh is the distance from the conical point along 
the /i-direction and for simplicity we have shown the 
dependence on only one of the two position variables 
(h,v). Here, the crystal shape has been oriented such 
that F(h, v) is constant along the first-order edge be- 
tween facets. The second term then corresponds to the 
jump in slope at the conical point. The third term shows, 
however, that there is more fundamental nonanalytic be- 
havior, namely that an analytic expansion about this 
point does not exist, even if one requires that the ex- 
pansion only be valid in the curved region of the ECS 
/ free energy surface. (Analogous nonanalyticities ap- 
pear in the surface tension F(x,y) p2].) Nonanalytic- 
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ities in the free energy [beyond the expected disconti- 
nuities in the derivatives of F(h, v)] at first-order coexis- 
tence boundaries have been discussed previously [|6| and 
are thought to be a rather generic feature in phase dia- 
grams. They are absent within mean field theory, where 
first-order transitions are associated simply with cross- 
ings of the local minima in the free energy surface. It 
is only through the effects of renormalization that such 
nonanalytic behavior can appear. In the Ising model, the 
nonanalyticities in the free energy F(H) at zero field (H 
= 0) for T < Tc are weak essential singularities. For 
vector spin models, as here, power-law singularities are 
observed. These divergences can be explained using spin 
wave theory (46| . 

Clearly, the question which must be asked is how gen- 
erally a conical point occurs, as opposed to having a ridge 
extend into the rounded region of the ECS. We are not 
able to answer this question rigorously. However, we be- 
lieve that the evidence from the six-vertex model sug- 
gests the conical point will be a general feature of equi- 
librium crystal shapes at the point where a sharp edge 
between facets meets a rounded region of the ECS. This 
belief is based upon the following observation: In the six- 
vertex model, the appearance of the conical point, hav- 
ing such a high degree of symmetry, does not seem to be 
due to any special symmetry which already exists in the 
original Hamiltonian. Rather, the symmetry seems to be 
generated spontaneously under renormalization (and this 
occurs for all temperatures and values of the interaction 
parameters that put us in the "low-temperature ferro- 
electric regime" of the model). Therefore, we expect such 
conical points to be a rather generic feature of equilib- 
rium crystal shapes, probably occurring for simple cubic 
as well as FCC models and outside the RSOS approxima- 
tion. On the other hand, Neergaard and Den Nijs have 
recently presented an argument that the conical point, or 
at least the specific scaling behavior at this point, may 
be a rather special feature of the six-vertex model p7| . 
In light of these conflicting opinions, we must say that 
the generality of this feature remains an open question. 

It is interesting to note that the conical point in the 
crystal shape and the associated coexistence region in the 
space of surface orientations is precisely the converse of a 
cusp in the surface tension and the associated facet on the 
crystal shape ]is| ]. In this latter case, the conical point 
(cusp) in F(x,y) corresponds to a flat facet in F(h,v). 
The analogue of the one-parameter family of orientations 
coexisting for the former case is the one-parameter fam- 
ily of points on the crystal shape (namely, those on the 
boundary of the facet) which all correspond to the cusp 
point in F(x,y). Also, the sharp edge, or ridge, between 
the facets in F(h, v) is analogous to the "groove" in the 
surface tension F(x,y) which exists between the two co- 
existence regions associated with the two conical points 
II- 

We are now ready to discuss the first question posed 
at the beginning of Section III B , namely the determina- 



tion of the boundary between stable and unstable orien- 



tations as a function of temperature. As noted above, 
this boundary can be determined exactly for the FCC 
case, within the [110]-RSOS approximation, because of 
the mapping onto the exactly-solved six-vertex model. 
For the simple cubic case, we do not have an exact so- 
lution even within the [lll]-RSOS approximation. How- 
ever, we can calculate the "opening angle" 9 C {T) for this 
curve. [In Fig. ^(b), this is the opening angle between the 
bold solid line and the octant boundary] This angle is 
also half the opening angle 28 C (T) between the two facet 
boundaries on the ECS at the conical point. Because 
calculation of 9 C (T) involves surfaces which are vicinal 
to the facet orientations (i.e., the steps across the facets 
are widely spaced so that step-step interactions are not 
important), it can be determined from an exact calcu- 
lation of the free energy of an isolated step within the 
[lll]-RSOS approximation. The calculation of the step 
free energy, f s (T, 9) and the angle C (T) are discussed in 
the appendix. An alternate way to view C (T) is that 
its inverse function, T ps (6), gives the "phase separation 
temperature" for a vicinal surface (i.e., that tempera- 
ture below which the surface is unstable) as a function 
of the angle of the steps across the surface (and thus its 
orientation). A plot of T ps (9) is given in Fig. |?]. 

We have carried out some simulations of the breakup of 
various surfaces [hkl] with h > k > I > at various tem- 
peratures. Those surfaces studied which are in the limit 
I < h ([15, 5, 1], [22, 7, 1], and [30, 29, 1] surfaces) appear 
to decompose into two phases, one vicinal to the [100] sur- 
face and the other vicinal to the [010] surface [see, e.g., 
Fig. ^(b)] . When the system is quenched from high tem- 
peratures, thin stripes of these two vicinal phases form 
and then coarsen over time until only one stripe of each 
phase remains in the system. We have measured the ap- 
proximate angles of the steps across these vicinal phases 
and find that they are in quite good agreement with our 
formula for 8 C (T) [see Fig. [|. 

Note that the fact that the surface appears to decom- 
pose into only two surface orientations, rather than an 
entire one-parameter family, seems to be in contradiction 
to the conical point scenario we have argued is likely to 
apply. However, three important points should be stated. 
The first is that it is, in fact, hard to tell if there are really 
only two orientations present. What we can say with con- 
fidence is that there seems to be a breakup into a "phase" 
with orientations near [100] and a "phase" with orienta- 
tions near [010], and that both these phases consist of 
vicinal surfaces with steps having approximately the an- 
gle C (T) across them. This is not really in contradiction 
with what we expect: In the conical point scenario, the 
decomposition will consist of a weighted distribution of 
the surfaces in the one-parameter family of marginally 
stable orientations. For decomposition of a surface close 
to being a [hkO] surface, the decomposition must primar- 
ily consist of those surfaces vicinal to the [100] and [010] 
surfaces in order to get the correct average surface orien- 
tation. (For a [hkO] surface itself, the decomposition will 
simply be into the [100] and [010] surface orientations.) 
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Our second point is that when we then carry out sim- 
ulations of surfaces which are not as close to being an 
[hkO] surface (and at temperatures nearer to Xcr in or- 
der that we see any phase separation at all), then the 
equilibrium configurations do not appear to be so simple. 
Fig. H shows an equilibrium configuration of a [13, 13,4] 
surface at T = 8J2. At T = 9J2, the configurations are 
even more complicated and it is in fact hard to tell if 
the surface has phase-separated at all. Unfortunately, in 
studying these configurations, we have found it difficult 
to quantify the surface phase separation since one would 
have to measure the surface orientation coarse-grained 
over some region, because (unlike, e.g., in the case of 
an XY model) the order parameter is not a continuous 
variable at the microscopic level: each rhombic tile rep- 
resents either a [100], [010], or [001] orientation. The 
correct scale over which to coarse-grain is not obvious, 
and of course, the larger this scale, the larger the system 
size should be in order to avoid finite-size constraints. 

Our final point is that the periodic boundary condi- 
tions put some constraints on the system, specifically, 
that the stripes of each "phase" in Fig. [|(b) must be 
horizontal. This may be influencing the results consid- 
erably. Therefore, one might need to go to much larger 
system sizes (which is not computationally feasible at 
present) in order to start truly seeing the other available 
surface orientations. Another approach would be to try 
to choose alternative boundary conditions (such as free 
or helical). However, implementing such boundary con- 
ditions and forming the initial surface configuration is a 
nontrivial task in this model. 

In the final analysis, we would have to say that the 
simulations currently do not provide conclusive evidence 
either for or against the conical point scenario in this 
model. 



2. Kinetics 

To study the kinetics of the phase ordering process, we 
turn to the continuum formalism used by Liu and Metiu 
[ p~5|] . For the case of the evaporation-condensation mech- 
anism, and neglecting noise since temperature is believed 
to be an irrelevant variable, they derived the Langevin 
equation 



dt 



-A 



aVV(r,t) - V(V 



dF {x,y) , 
dip ' 



, (8) 



where r = (h, v) and ip = (x,y). (Here, we have used 
the constraint V x ip = to rewrite the first term on the 
right hand side ||] . This constraint follows from the fact 
that xjj is the gradient of the surface.) If surface diffusion 
is the dominant mechanism, then there is an additional 
— V 2 in front on the right hand side. Fo(x,y) is some 
suitably coarse-grained projected surface tension, not the 
full thermodynamic projected surface tension F(x, y). In 
particular, Fq(x, y) is not necessarily convex since, in the 



case when the (x, y) orientation is unstable, it reflects the 
free energy of metastable configurations having orienta- 
tion (a;, y) rather than the thermodynamic configurations 
once phase separation has occurred (See, e.g., Ref. p9|). 

Because the order parameter ip is conserved, the stable 
stationary states {i.e., those orientations into which a 
surface decomposes) are not simply given by the minima 
of Fo(x,y). Rather they are given by those values of i/> 
that share a common support tangent on Fo(x,y) (See, 
e.g., Fig. 3 in Ref. ||). 

Note that the structure of Eq. (|^) is very similar to 
that for the traditional phase separation process, i.e., 
the so-called Cahn-Hilliard Equation (See Q). The one 
difference is the reversal of the ordering of the two vec- 
tor operators in the second term on the right hand side: 
Here we have the gradient of a divergence (and also the 
constraint that V X ip = 0), while in the traditional case 
one has the divergence of a gradient Q . The former case 
causes a coupling of the different components of the or- 
der parameter, and of the spatial and order-parameter 
degrees of freedom, while the latter case leaves them un- 
coupled. (In fact, the former case only makes sense if 
n = d where d is the spatial dimension of the interface 
and n is the number of degrees of freedom of the order 
parameter.) 

Before considering the case of most interest to us at 
present (where phase separation is expected to occur into 
a continuous family of orientations), we will first discuss 
the "simpler" case of separation into a finite number of 
orientations. This has been the subject of much very re- 
cent work. In Ref. | |L5| , the authors assumed a form for 
the surface tension which results in a few (e.g., three) 
stable stationary states in Fq(x, y). On the basis of simu- 
lation results and power-counting, they argued that the 
growth of the domain size likely satisfies L(t) ~ i 1 / 4 . 
However, in the case of traditional phase separation (i.e., 
governed by the Cahn-Hilliard Equation) into a finite 
number of phases, such power-counting is known to fail. 
Because of the sharpness of domain walls, one finds in- 
stead L(t) ~ t 1 / 3 . The most general arguments for these 
results are given by Bray Q , who argues that the growth 
law in the case of conserved dynamics can be determined 
from the energy cost of a domain wall. If the energy cost 
goes like E(L) ~ V then L(t) ~ *V(2+<*-iO. For the case 
of sharp domain walls, e.g., in an Ising model, y = d— 1; 
whereas in the case of domain walls which themselves are 
spread out over a width L, as in the isotropic vector spin 
models, y = d — 2. 

Do Bray's arguments continue to hold for Eq. (||), or 
does the reversal in the order of the vector operators re- 
sult in different behavior? It appears that, in fact, the 
renormalization flow equations change such that simple 
connection between y and the growth law exponent no 
longer necessarily follows j|. Indeed, in addition to Liu 
and Metiu |L5| |, a few other groups have recently sim- 
ulated equations similar to Eq. (0) and have found re- 
sults most compatible with L(t) ~ t 1 ^ 4 ||||. Thus it 
seems possible that the order of the vector operators does 
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change the kinetics of the phase separation process. On 
the other hand, given the history of low measured expo- 
nents (often around 1/4) in Ising systems where eventu- 
ally the exponent was shown to cross over to 1/3 [pl| , 
we believe that, in the absence of further analytical un- 
derstanding, both L(t) ~ t 1/4 and L(t) ~ t 1 / 3 (or even 
some intermediate exponent!) should still be considered 
as viable candidates for the asymptotic behavior. 

Finally, we note that we have also published some of 
our own numerical evidence bearing on this question. In 
Rcf . [|l3| , we looked at the behavior of the characteristic 
length scale for the coarsening of a [111] surface cooled 
slowly at a constant rate T. Assuming the "underlying 
growth law" in the absence of diverging barriers would be 
L(t) ~ i 1 / 3 , we argued that the behavior under asymp- 
totically slow cooling would be L(T,T = 0) ~ F' 1 / 4 . 
It was noted that this expectation was in reasonable 
agreement with our numerical results. Here, we add two 
comments: (1) If the "underlying growth law" is instead 
L(t) ~ t 1 / 4 , then the corresponding behavior at a con- 
stant cooling rate is L(T,T = 0) ~ r -1 / 5 (See endnote 
68 of Rcf. |L3|]). (2) The simulations performed there, 
and further simulations carried out on systems as large 
as 240 2 p3], cannot conclusively distinguish these two 
possibilities, with the best estimate of the exponent of T 
being —0.23, i.e., about halfway between the two expec- 
tations. 

Of course, the considerations of the previous three 
paragraphs apply to the case where the surface tension is 
such that phase separation occurs into a finite number of 
different orientations separated by sharp domain walls. 
This will be relevant to our discussions in Sec. IV. In the 



present case, however, our Fo(x,y) is expected to be of 
such a form that we have phase ordering into a continu- 
ous one-parameter family of orientations. The numerical 
work of Siegert and Plischke || suggests that in this case, 
as in the standard case when the vector operators are re- 
versed (e.g., an isotropic XY or Heisenberg model), the 
exponent for the growth law is 1/4. (Note however that, 
technically speaking, the case n = d = 2 is excluded 
from Bray and Rutenberg's theory of phase ordering in 
the Cahn-Hilliard equation because of l ong -range corre- 
lations between the topological defects |5C]] . Simulations 
of this case are nonetheless consistent with L(t) ~ t 1 ! 4 , 
but perhaps with logarithmic corrections and deviations 
from scaling p3fl.) 

It is worth emphasizing that the whole decomposi- 
tion process here is quite different from the case where 
separation occurs into a finite number of orientations: 
Here, there are really no well-defined domains or "facets" 
formed at all; rather, it is the length scale over which the 
surface bends that will increase over time (which can be 
measured most easily from a correlation function or its 
Fourier transform). Thus terms such as "faceting" and 
"phase separation" seem somewhat inappropriate. 



IV. CRYSTAL GROWTH OR ETCHING 

We now consider the case in which the crystal is 
slightly out of equilibrium with its vapor, that is, the case 
of slow growth or etching. In fact, it is known that, far 
out of equilibrium, the nonequilibrium effects can drive 
the formation of facets even when such faceting would not 
occur in equilibrium |5]](|. Here, however, we will assume 
that we are close enough to equilibrium that the stability 
of a surface is still determined from equilibrium consid- 
erations, and will study the effect of a chemical potential 
difference A/i between solid and vapor phases on the 
dynamics of the phase separation. In the magnetic lan- 
guage of the 3-D Ising ferromagnet, such a situation can 
be thought of as the application of a uniform magnetic 
field in the model. 

First we consider the temperature regime T < Tqr 
and, for definiteness, we discuss the case of crystal growth 
(A/i < 0). Then the free energy cost associated with 
adding surface atoms to produce a step across a facet 
can be written approximately as follows: 



AF^f p (T,6)L 



L 2 

|A/i| — tan(0) 



(9) 



Here f p (T,6) and L are the projected step free energy 
and the projected length of the step, respectively (see 
appendix). The projection is onto one of the lattice axes, 
and 9 is the angle of the step relative to that axis. [Eq. 
(^) is valid when the linear size of the facet itself is larger 
than L.\ The structure of the problem is the same as 
that of determining the growth rate (due to nucleation) 
of a growing crystal surface at a temperature below its 
roughening transition j36|: Assuming that |A/i| is not too 
large, then for small L the first term in the expression 
dominates and AF increases linearly with L. However, 
for large enough L, the second term dominates and AF 
is a decreasing function of L. 

To find the smallest possible barrier associated with 
propagating a step across an infinite facet, we want to 
find the length of step Lq and the angle 9q such that 
the free energy cost is maximized as a function of L 
and minimized as a function if 6. This we do by setting 
d{AF) /dL = d(AF) /d9 = Q. We find 



^ f P (T,0 o ) 
~ |A/i|tan(0 o ) ' 

with an implicit expression for 9q given by 



df P (T, 



d[tan(0)] 



fp(T,9 ) 
2tan(6> ) 



(10) 



(11) 



However, this latter expression always has the solution 
tan(0o) = 1, independent of A/i or T, so a step of 45° 
always gives the smallest barrier. Thus, ( JTTj| ) becomes 



L nf p (T,e = ir/4)/\An\ 



(12) 







and the free energy barrier is 



/ p 2 (r,fl = 7 r/4) 

2|A M | 



(13) 



Eq. (|T^) gives the maximum size of the barriers for a 
coarsening surface. The barriers to coarsening will grow, 
and thus the coarsening of the surface will look roughly 
logarithmic, out to a time of order 



exp[L 2 (T,0 = 7r/4)/(2|A M |T)] 



(14) 



when the characteristic facet size L is of order Lq. At 
times longer than this, the barrier height to coarsening 
remains saturated at AF max , independent of L, and the 
growth will thus proceed with a power law 



L(t) ~ (7t)» 



(15) 



with 7 w exp[-/|(T,0 = 7r/4)/(2|A/i|T)]. Here, n is 
the exponent for growth associated with Eq. (|§|) in the 
case where the surface breaks up into a finite number of 
orientations (thus, either n = 1/3 or 1/4, as discussed in 
Section p!B2| ). 

The form of 7 is closely analogous to that for the rate of 
growth of a defect-free crystal surface below the roughen- 
ing transition p6| . The one important difference, how- 
ever, is that in this latter case, no driving force exists 
for the growth of the crystal in the limit A/i = 0. In 
the case of the phase separation problem, however, while 
there is still no driving force for growth of the crystal in 
the absence of a chemical potential difference, there is a 
driving force for phase separation of the surface. Thus, 
the characteristic length scale grows with time even when 
A/i = 0, but the process occurs logarithmically-slowly in 
this case. 

We have performed Monte Carlo simulations of the 
[111J-RSOS model in order to test certain aspects of our 
above arguments. First, in order to study the tenet of 
our argument involving the behavior of the free energy 
barrier as a function of the size of the facet [Eq. (g)], we 
consider the process of removing one face of a "cubical 
projection" from a crystal surface, as illustrated in Fig. 
^. Fig. |l^ is an Arrhenius plot of the time such a process 
takes as a function of the linear size L of the projection 
for the case Afj, = 0. (Cf. Fig. 6 of [|l3| where a similar 
plot was presented for the case of fully three-dimensional 
coarsening. Also see that reference for further details 
about the simulation methods.) We see that, as is our 
expectation, the slope on the Arrhenius plot, and thus 
the free energy barrier, is an increasing function of the 
size of the projection. In fact, the expected free energy 
barrier in the limit T — ► can easily be calculated |5J] 
to be 



F B (L,T = 0) 



12 J 2 , L = 2 

4J 2 (L + 2) , L > 2 



(16) 



Fits to the form t = t (L)e FB{ - L ' T = 0)/T , with t (L) as 
a free parameter, are also shown in the figure and are 



generally very good for the low-temperature data. Only 
for L = 8 are systematic deviations from the fits still 
evident at the lowest temperatures, with the slope on the 
plot being greater than that predicted. This deviation is 
due to our approximation of Fb(L,T) by Fb(L,T = 0) 
(See Sec. II of Ref. f|). 

In Fig. [H], we present an Arrhenius plot for the time 
to shrink a cubical projection of size L for A/i/J 2 = 2. 
(Here, the sign of A/x is chosen to be positive, i.e., to 
favor shrinking of the projection). The contrast to Fig. 
|l0| is apparent: For small L, the free energy barrier is an 
increasing function of L. However, for L > 6, the barrier 
appears to saturate. The expected barrier for each L can, 
in fact, be calculated by considering the highest energy 
state the system must pass through in removing a face 
on the projection. This gives 



F B (L,T = 0) = 



12J 2 - A^ , L = 2 . 

4J 2 (L + 2) - Afi(L + 1) , 2 < L < 6 . 
18J 2 , L > 6 . 

(17) 



(Here, the saturation of the barrier at L — 6 is for 
the specific case of A/i/J 2 = 2.) Fits to the form 
t = < (i)e Fs(i ' T=0)/T , with t (L) as a free parameter, 
are also shown in the figure and are generally very good 
for the low-temperature data. (The one fit which shows 
a bit of deviation is for L — 6, probably because of the 
presence of a large degeneracy of barrier states with very 
similar barrier heights.) Also, the value of 18J 2 which 
we get for the saturated barrier height is in quite good 
agreement with the estimate of 16 J 2 provided by Eq. Jl3| ) 
in the limit T — > [which can be worked out using Eq. 
( |i"l| ) and formulas in the appendix]. 

The preceding study tested our assertions about the 
barrier heights to coarsening for one very artificial kind 
of configuration. In order to test more directly the effect 
of a nonzero A/i on the coarsening process, we present 
in Fig. [l^ full-scale simulations of the coarsening of a 
[111] surface, which is initially in a random (T = 00) 
configuration and is quenched to T — 2J 2 . The figure 
shows the growth of the characteristic length scale L(t) 
over time for various values of | A/x|/ J 2 . We see that, in its 
qualitative aspects, the Monte Carlo data is in agreement 
with our expectations: First, we consider the case A/j, = 
0. The growth, although not yet quite logarithmic as 
would be expected at asymptotically long times, is slower 
than power law; the effective exponent in L(t) ~ i" off is 
n e s ~ 0.07 at the latest times and is continuing to slowly 
decrease. 

As the chemical potential difference between solid and 
vapor is increased from zero, the growth becomes faster. 
For I Afj,\/ J 2 — 0.5, the effective exponent is virtually con- 
stant over the six decades of time shown. For the larger 
\Afi\, the increase in the effective exponent with time 
is quite gradual. However, the most marked bending of 
the data does appear to occur at times whose trend with 
|A/i| is in qualitative agreement with (but that seem to 
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be quantitatively a bit later than) the expectations of Eq. 
(fj). For |A/x|/J 2 = 1.2 and 2.0, there is a regime where 
the effective power law approaches (although remains a 
bit less than) 1/4. 

Also visible for A//|/ J2 = 2.0, and just barely for 
| A/i|/ J2 = 1.2, is a rather sharp leveling off in the growth 
of L{t) at late times. Studying the configurations at 
these late times, we find that the system appears to have 
reached a steady-state in which at any time, there are 
structures such as ledges, etc. because the crystal is in 
the process of growing (or shrinking), due to the nonzero 
value of A/i. The fact that such growth is in progress puts 
a maximum limit on the average facet size. [Our measure 
of L(t) as being inversely proportional to the total length 
of boundary between the different facet orientations may 
also be underestimating somewhat the true value of the 
characteristic length scale. See Sec. IV. C of Ref. [fl3f .] 
A more detailed study of correlation functions might be 
useful in further characterizing this behavior (e.g., one 
might even find a breakdown of scaling). Clearly, these 
are effects which occur in a regime beyond which our 
assumption of being very near equilibrium is valid. 

In closing this section, we briefly discuss the likely ef- 
fects of a nonzero chemical potential difference in the 
temperature regime Tcr < T < Ter. Here, the coarsen- 
ing is already expected to be a power law for A/j, = 0. A 
nonzero A/x seems unlikely to do more than alter the 
prefactor of the power law. One important question, 
however, is whether it might provide a relevant perturba- 
tion which will break the symmetry associated with the 
conical point and thus allow the selection of two partic- 
ular orientations for the surface to separate into. In the 
absence of a renormalization group picture for the coni- 
cal point on the ECS, we do not know how to determine 
if this perturbation is relevant, and hence this question 
remains unanswered. 



V. SUMMARY AND CONCLUSIONS 

In this paper, we have considered the "phase separa- 
tion" of crystal surfaces within an Ising lattice gas model 
for materials with the sodium chloride structure. We 
found that, depending on the temperature and growth 
conditions, a number of interesting behaviors can be ob- 
served. 

For a crystal in equilibrium with its vapor, we argued 
that the coarsening of the surface structure will be loga- 
rithmic in time for temperatures below the corner round- 
ing transition temperature Tcr, where an arbitrary sur- 
face [hkl] (with h, k, I > 0) undergoes three-phase sepa- 
ration into the smooth [100], [010], and [001] facets. Such 
sluggish dynamics occur because the coarsening process 
involves the creation of steps across the smooth facets 
and these steps have a nonzero free energy per unit length 
below Tcr. Since the discreteness of the system plays 
a fundamental role here, recently proposed continuum 



models for the dynamics of phase separation of surfaces 
fail to capture this effect. 

At all temperatures between Tcr & n d the edge round- 
ing temperature Ter, some surface orientations remain 
unstable. It appears likely, based on the exact solution 
of a related model, that (because of a spontaneously- 
generated symmetry) such surfaces will generically un- 
dergo decomposition into an entire one-parameter fam- 
ily of surfaces. This process, analogous to phase or- 
dering in an XY (or Heiscnbcrg) model, should obey 
L(t) ~ t 1 / 4 . Only if effects not considered (such as elastic 
interactions or surface melting) break the spontaneously- 
generated symmetry, does the possibility exist for simple 
two-phase separation. In that case, the growth law will 
be L(t) ~ t n , with either n = 1/3 or 1/4. 

Such an algebraic growth law should also hold at late 
enough times below Tcr for the case of driven surfaces, 
i.e., where the crystal and vapor are not in equilib- 
rium (crystal growth, or etching), with the crossover time 
from very slow (approximately logarithmic) to power-law 
growth diverging rapidly as the chemical potential differ- 
ence between crystal and vapor decreases to zero. We 
noted the close analogy in this case to the growth of a 
crystal facet by nucleation below the roughening transi- 
tion. 

In closing, we should mention some effects which have 
been neglected in our microscopic model and could con- 
ceivably be important in real materials. First, there 
are impurities and defects which may act to speed up 
the phase separation of surfaces, just as they do for the 
growth of crystal surfaces below the roughening transi- 
tion [ (36|]55f| . This could be particularly important given 
that impurities often tend to segregate to the surface, so 
that even reasonably pure samples can have fairly dirty 
surfaces. Second, there are also elastic effects, which have 
been discussed in recent work JlJJl^,^,^] as possibly 
slowing down the faceting process and even stabilizing 
the system at a maximum facet size. There is some re- 
cent experimental evidence which seems to support this 
scenario [[l0[. Finally, there is the possibility of surface 
melting |38| , [48| and surface reconstructions |L2| which can 
result in a more complicated ECS. Clearly these issues 
must all be addressed in order to obtain a realistic pic- 
ture of the thermal faceting in real materials. 
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APPENDIX: CALCULATION OF THE STEP 
FREE ENERGY AND RELATED QUANTITIES 

Here, we calculate the step free energy, f s (T, 9) , as a 
function of temperature T and step angle 9 (measured 
with respect to one of the lattice axes) within the [111]— 
RSOS approximation. We then use this result to discuss 
features of the phase separation problem in the temper- 
ature regimes both below and above Tcr. The reader is 
also referred to Ref. where a calculation of the step 
free energy was performed from a somewhat different, 
but related, perspective. 

It is in general quite difficult to compute the step free 
energy in a so-called "canonical ensemble," where we re- 
strict the step to have a certain angle. The easier route is 
to compute the step free energy in the "grand canonical 
ensemble" , where we do not constrain the step angle, but 
instead apply a "field" h which couples to the angle of 
the step and thus allows us to control the average step 
angle |5q] . To obtain f s (T, 0) from f s (T, h) we then need 
only determine this average step angle 9 as a function 
of h. Then, f s (T, 9) is related to f s (T, h) by a Legendre 
transform, which has an intuitive physical interpretation. 

In the [lll]-RSOS approximation, a step across, say, a 
[100] surface consists simply of a one-dimensional series 
of plaquettes of either the [010] or [001] orientation [see, 
e.g., Fig. ^(b)]. If we apply a "field" h in order to favor 
one orientation over another, then the Hamiltonian for 
the step is equivalent to that for a one-dimensional Ising 
model in a magnetic field. Within this grand canonical 
ensemble, there are several routes to calculate the step 
free energy. If one makes the analogy to the Ising model 
then one can write down the free energy per plaquette 
[ 20 1 . This is done most easily by using a transfer matrix 
[ 59 1 . However, to make connections to free energy barri- 
ers, we find it most natural to compute the free energy 
per unit projected length along one of the lattice axes. 

The calculation proceeds as in Sec. II of Ref. jl3| , ex- 
cept with the addition of a field term —niih to the energy 
of a column in Eq. (2.10) of p3[ . Here, we are interested 
only in the large-size (L — > oo) limit in which we can 
write the partition function as 



Z{L) = 



_ p -4J 2 L/T 



1 



e -(8J 2 -h)/T 
1 _ e -(4J 2 -?0/T 



(Al) 



where L is the projected length along the [010] axis. The 
quantity f p (T, h), a free energy per unit projected length 
at a given applied "field" h is then 



f p (T,h)= -jln[Z(L)] 
= 4J, - Tin I 1 



q 2 H 
1-qH 



(A2) 



where q = exp(— 4J 2 /T) and H = exp(h/T). The (aver- 
age) tangent of the angle of the step with respect to the 
[010] axis is given by 



tan(0)= - 



df P (T, h) 
dh 



(1 - qH)(l - qH + q 2 H) 



(A3) 



This equation can be inverted in order to express h as a 
function of r = tan(#): 



Tin 



(2 - q)r + q - y/(r 2 + l)q 2 + 2(2 - q)qr 
2q(l-q)r 



by 



(A4) 

The step free energy as a function of angle is then given 



fs(T, 6) = [f p (T, h) + h tan(0)] cos(#) , (A5) 



where, for a given angle 9, h is determined by (A4). 
The addition of hta,n(9) takes the Legendre transform 
of f p (T,h) to give us f p (T,9). The physical interpreta- 
tion of this transform is simply that this term subtracts 
off, from f p (T,h), the energy contribution due to the 
coupling of the field to the angle, thus leaving us with 
what the free energy of a step at angle 9 would be in the 
absence of any applied "field." The factor cos(#) is neces- 
sary to give us the free energy per unit (macroscopic) step 
length, f s , rather than the free energy per unit projected 
length, f p . 

There are two ways in which the calculation above is 
useful to us. First, for T < Tcr and Afi — 0, the free 
energy barrier per unit length for shrinking of a cubi- 
cal projection of edge length L (in the large L limit) is 
given rigorously by fs{T) = f P {T, h = 0). This is also, 
roughly speaking, the appropriate barrie r to c onsider for 
a coarsening surface (as discussed in Sec. Ill A), although 
the value of /s(T") can no longer be so rigorously deter- 
mined. In particular, it would depend on precisely how 
we choose to define the characteristic length scale L. For 
coarsening in the case Afj, ^ 0, the important quantity 
that enters into the equations is f p (T,9 — 7r/4) as dis- 



cussed in Sec. IV 



Second, we can extend the above results to derive the 



"opening angle" 9 C (T), discussed in Sec. Ill B, in the tem- 
perature regime Tcr < T < Ter. To determine this 
angle, we consider a surface [hkl] with < k < h and 
< I <C h. Such a surface is thus very close to being 
a [hkO] surface and if it is unstable then (independent 
of whic h of th e conical point or ridge scenarios discussed 
in Sec. 



[II B 



is correct), it will break up (at least pri- 
marily) into surfaces which are vicinal to the [100] ori- 
entation and surfaces which are vicinal to the [010] ori- 
entation. These surfaces will have steps across them of 
a definite angle, given by 9 C (T). This opening angle is 
then determined by the condition that it minimize the 
total free energy of the resulting surface subject to the 
constraint that the surface has the correct average tilt in 
the z direction. The free energy of the surface is given 
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by f s (T,9)L s where L s is the total length of step. The 
tilt is given by L s sm(0). Therefore, C (T) is that angle 9 
which minimizes the quantity f a (T, 6)/ sin(0) JS0| . This 
one-dimensional minimization can easily be performed 
numerically. The results are shown in Fig. ffl. 
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FIG. 1. The decomposition of a surface in the [111J-RSOS 
model for the Hamiltonian of Eq. ([j]) at T < Ten,. Shown is 
a [111] surface which has been quenched from infinite tem- 
perature to T = 3J2 at times (a) t = 0, (b) t — 100, and (c) 
t = 10000 (in MC steps/plaquette) following the quench. 

FIG. 2. Qualitative thermal evolution of the equilibrium 
crystal shape (ECS) for the model specified by the Hamilto- 
nian of Eq. @), after Ref. [jl|. 

FIG. 3. The stable surface orientations for the model spec- 
ified by the Hamiltonian of Eq. ([j]) in various temperature 
regimes. This figure shows an eighth of a sphere, represent- 
ing surface orientations by the unit vector normal to the sur- 
face. The shaded regions on the figure show those surface 
orientations which appear on the ECS of Fig. |^. All surface 
orientations in the unshaded regions are thermodynamically 
unstable |tJ and surfaces prepared in such orientations will 
decompose into some combination of stable surfaces. 
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FIG. 4. Equilibrium configurations for a [15,5,1] surface in 
the [lll]-RSOS model at various temperatures: (a) T = U2, 
(b) T = 14 J 2 , and (c) T = 34 J 2 . Minimization of the total 
surface free energy, under the assumption that the surface is 
close enough to being vicinal that the steps do not interact, 
gives r ps ([15, 5, 1]) w 22.6J2 for the phase separation temper- 
ature (See appendix and Fig. |^) . This estimate should provide 
an upper bound for the actual value of T ps ([15, 5, 1]). 

FIG. 5. The "ridge scenario." In (a), we show a closeup 
of the corner region of the ECS for Ten, < T < Ter, un- 
der the assumption that the sharp edges between the facets 
continue into the curved part of the ECS, finally terminat- 
ing at second-order critical points. At any point along this 
edge, there is "coexistence" between two surface orientations. 
(Here, sharp first-order edges are shown by dashed lines while 
the second-order Pokrovsky-Talapov boundaries between the 
facets and the curved part of the ECS are shown by solid lines. 
The solid circles represent second-order critical points at the 
end of the first-order lines.) In (b), we present the diagram 
of the stable surface orientations, showing how any unstable 
surface will phase separate into the two coexisting stable sur- 
faces. In our example, any orientation along the dashed curve 
(such as the one indicated by the open circle) decomposes into 
the two (marginally) stable surface orientations indicated by 
the solid circles. 

FIG. 6. The "conical point scenario." In (a), we show the 
corner region of the ECS for Tor < T < Ter, under the as- 
sumption that the sharp edge between the facets terminates 
at the point where it meets the curved region of the ECS. 
At this "conical point" (indicated by a solid circle), an entire 
one-parameter family of surface orientations coexist. In the 
curved part of the crystal shape, these orientations come to- 
gether at this point just as they do at the tip of a cone. In 
(b), we present the diagram of the stable surface orientations, 
showing how an arbitrary unstable surface orientation (shown 
by an open circle) will decompose into a combination of the 
entire one-parameter family of coexisting surface orientations 
(shown by the bold solid curve). 

FIG. 7. Phase separation temperature T ps (0) for a vicinal 
surface as a function of the angle of the steps across it in 
the [lll]-RSOS model, as calculated in the appendix. Alter- 
natively, as explained in the text, the inverse function 8 C {T) 
can be thought to represent the angle of the steps across the 
vicinal surfaces which an unstable crystal surface [hkl] with 
< k < h and < I <C h decomposes into. The symbols 
represent measurements of C (T) from simulations of this de- 
composition process. For a finite value of J1/J2 (i.e., relaxing 
the RSOS constraint), the curve would look qualitatively sim- 
ilar for large angles but would not diverge as 6 —* 0. Instead, 
we would have T ps (6 — » 0) = Ter. 



FIG. 8. An equilibrium configuration for a [13,13,4] surface 
in the [lll]-RSOS model at T = 8J 2 . Note that the surface 
appears to have phase-separated, but the equilibrium state is 
quite complicated. This may be an indication that many dif- 
ferent surface orientations are involved, as one would expect 
from the conical point scenario. 

FIG. 9. Shrinking a cubical projection in the [lll]-RSOS 
model, (a) shows a cubical projection with L = 4. (b) shows 
the same projection once one face has been removed. For 
A/i > 0, this configuration is lower in energy than that of (a) 
by 4J 2 + A/iL 2 . 

FIG. 10. The time to shrink a cubical projection for 
Afj, = 0. Shown are the results of Monte Carlo simulations for 
the average time to remove the first face from a cubical pro- 
jection of size L [such as that shown in Fig. |^(a)] for the case 
Afj, = 0. Each point is an average over 900 runs with stan- 
dard error smaller than the symbol size. The dotted lines are 
one-parameter fits to the form t = to(L)e FB(L ' T ^°^ T with 
Fb(L,T — 0) given by Eq. ( |l6| ) and to(L) a free parameter. 
Between 4 and 15 of the lowest temperature data points are 
used for each fit. 

FIG. 11. The time to shrink a cubical projection for 
nonzero Afi. Shown are the results of Monte Carlo simu- 
lations for the average time to remove the first face from a 
cubical projection of size L for the case A/i/J 2 = 2. Each 
point is an average over 900 runs with standard error smaller 
than the symbol size. The dotted lines are one-parameter fits 
to the form t = t Q e FB{L ' T = 0),T with F B (L,T = 0) given by 
Eq. (Jl?]) and ta(L) a free parameter. Between 6 and 18 of 
the lowest temperature data points are used for each fit. For 
L < 4, some of the fitted data lies off scale. 

FIG. 12. Growth of the characteristic length scale L(t) dur- 
ing coarsening of a [111] surface at T = 2.J2 for various val- 
ues of Afj (given on the figure in units of J2). The runs 
were performed on a 120 2 system. For each value of Afj, the 
Monte Carlo data has been averaged over 77-150 runs, with 
error bars showing the standard error. A line of slope 1/4 
has been shown to facilitate comparison with the power law 
L(t) ~ t 1 / 4 . The characteristic time when we would expect 
a change in growth law on the basis of Eq. ( Jl4| ) and using 
f p (T = 2J 2 ,6 = vr/4) « 6.75J2, is to « 3 x 10 2 , 1 x 10 4 , 
2 x 10 6 , and 8 x 10 9 for A/x = 2.0, 1.2, 0.8, and 0.5, respec- 
tively. 
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